﻿#include <gtest/gtest.h>

#include "algo/algo_interpolate.h"
#include "algo/algo_projection.h"
#include "geom/geom_curve_bezier.h"
#include "geom/geom_curve_bspline.h"
#include "kernel/spline.h"

#include <algorithm>
#include <iostream>
#include <random>
#include <time.h>

#include <dbg.h>

using namespace CAGD;

TEST(TestInterpolate, BSplineCurve)
{
    std::mt19937 generator((unsigned int)time(nullptr));      // 使用 time(nullptr) 作为随机数生成器的种子
    std::uniform_real_distribution<double> uniform(0, 100.0); // 在 [0,1] 区间内的均匀分布

    std::vector<BPnt> points;

    int N = 100;
    for (int i = 0; i < N; i++)
        points.push_back({uniform(generator), uniform(generator), uniform(generator)});

    auto curve = Interpolate::BSplineCurve3(points);
    for (int i = 0; i < N; i++)
    {
        double t = 0.5;
        BPnt Q;
        auto distance = Projection::Near(curve, points[i], t, Q, CAGD_ABS_ACC);
        EXPECT_TRUE(distance < CAGD_ABS_ACC);
    }
}

TEST(TestInterpolate, BSplineSurface)
{
    auto f = [](float u, float v) -> BPnt { return BPnt{u, sin(u + v) * cos(u - v), v}; };

    // 20 x 30 个采样点
    int N = 10;
    int M = 15;
    std::vector<std::vector<BPnt>> points(N, std::vector<BPnt>(M));
    const double PI = 3.1415926535;
    for (int i = 0; i < N; i++)
        for (int j = 0; j < M; j++)
            points[i][j] = f(10.0 / (N - 1) * i, 10.0 / (M - 1) * j);

    // TIMER_START(InterpSurface)

    auto surface = Interpolate::BSplineSurface3x3(points);

    // TIMER_END(InterpSurface)

    double U1, U2, V1, V2;
    surface->Bounds(U1, U2, V1, V2);

    for (int i = 0; i < N; i++)
    {
        for (int j = 0; j < M; j++)
        {
            double u, v;
            u = (U1 + U2) / 2;
            v = (V1 + V2) / 2;
            BPnt Q;
            auto angle = Projection::Ortho2(surface, points[i][j], u, v, Q, CAGD_ABS_ACC);
            auto distance = Q.Distance(points[i][j]);

            // 最近距离投影不稳定
            // auto distance = Projection::Near(surface, points[i][j], u, v, Q, CAGD_ABS_ACC);
            // auto angle = (surface->Normal(u, v).Normalized() % (points[i][j] - Q).Normalized()).Modulus();

            // 考虑几种情况：角度正交，距离足够近；投影参数在边界上（每次迭代都出界的情况）
            bool pass = angle < CAGD_ABS_ACC || distance < CAGD_ABS_ACC || (fabs(u - U1) < CAGD_ABS_ACC * 2) ||
                        (fabs(u - U2) < CAGD_ABS_ACC * 2) || (fabs(v - V1) < CAGD_ABS_ACC * 2) ||
                        (fabs(v - V2) < CAGD_ABS_ACC * 2);
            EXPECT_TRUE(pass);
        }
    }
}